#We start by loading our packages

library(rayshader) 
library(rgl)
library(av)
library(gifski)
library(magick) #The first five are just for the 3d-model gif. 
library(ggplot2) 
library(dplyr)
library(magrittr) #should be included in dplyr, but didn't work?
library(plotly) #For interactive ggplot
library(hrbrthemes) #For better ggplot-themes
#First we load our data into our environment

library(tidytuesdayR) #Specific package for loading tidy-tuesday data

tuesdata <- tidytuesdayR::tt_load('2025-06-24') 

cases_month <- tuesdata$cases_month #Loading datasets into environment
cases_year <- tuesdata$cases_year   #Loading datasets into environment
#Looking at variable names
names(cases_month)
##  [1] "region"                "country"               "iso3"                 
##  [4] "year"                  "month"                 "measles_suspect"      
##  [7] "measles_clinical"      "measles_epi_linked"    "measles_lab_confirmed"
## [10] "measles_total"         "rubella_clinical"      "rubella_epi_linked"   
## [13] "rubella_lab_confirmed" "rubella_total"         "discarded"
names(cases_year)
##  [1] "region"                                                         
##  [2] "country"                                                        
##  [3] "iso3"                                                           
##  [4] "year"                                                           
##  [5] "total_population"                                               
##  [6] "annualized_population_most_recent_year_only"                    
##  [7] "total_suspected_measles_rubella_cases"                          
##  [8] "measles_total"                                                  
##  [9] "measles_lab_confirmed"                                          
## [10] "measles_epi_linked"                                             
## [11] "measles_clinical"                                               
## [12] "measles_incidence_rate_per_1000000_total_population"            
## [13] "rubella_total"                                                  
## [14] "rubella_lab_confirmed"                                          
## [15] "rubella_epi_linked"                                             
## [16] "rubella_clinical"                                               
## [17] "rubella_incidence_rate_per_1000000_total_population"            
## [18] "discarded_cases"                                                
## [19] "discarded_non_measles_rubella_cases_per_100000_total_population"
#Recoding regions for OUR understanding
cases_year <- cases_year %>% 
  mutate(region = recode(region, `AFRO` = "African Region", `AMRO` = "Regon of Americas", `EMRO` = "Eastern Mediterraenean Region", `EURO` = "European Region", `SEARO` = "Sout-East Asian Region", `WPRO` = "Western Pacific Region"))
cases_month <- cases_month %>% 
  mutate(region = recode(region, `AFR` = "African Region", `AMR` = "Regon of Americas", `EMR` = "Eastern Mediterraenean Region", `EUR` = "European Region", `SEAR` = "Sout-East Asian Region", `WPR` = "Western Pacific Region"))
unique(cases_month$region)
## [1] "African Region"                "Regon of Americas"            
## [3] "Eastern Mediterraenean Region" "European Region"              
## [5] "Sout-East Asian Region"        "Western Pacific Region"
#First we want to look at the overall average incidence rate pr. year worlwide
cases_year %>%
  group_by(year) %>% #Grouping by year
  mutate(avg_inc = mean(measles_incidence_rate_per_1000000_total_population)) %>% #Creating new variable with mutate
ggplot(aes(x = year, y = avg_inc)) + #inserting mapping-argument in overall ggplot for less typing
  geom_point() +
  theme_bw() +
  labs(title = "Avg. measles incidence rate worldwide \n(pr. 1,000,000)", x = "Year", y = "Avg. Incidence Rate") + #Notice the \n for changing                                                                                                                      #line
  geom_smooth()

#Now we want to facet by region. We supress warnings and messages in the output in the first row. 
plotfacet <- cases_year %>%
  group_by(year, region) %>%    #We group by year AND region so we get an annual avg. pr. region
  mutate(avg_reg = mean(measles_incidence_rate_per_1000000_total_population)) %>% 
ggplot() + 
  geom_point(aes(x = year, y = avg_reg, colour = region, shape = region)) +
  theme_bw() +
  labs(title = "Avg. measles incidence rate pr. region pr. year(pr. 1,000,000)", x = "Year", y = "Avg. Incidence Rate", color = "Region", shape = "Region") +
  geom_smooth(aes(x = year, y = avg_reg), color = "grey", se=FALSE, show.legend = FALSE) + 
  facet_wrap(~ region, scales = "free") +
  geom_line(aes(x = year, y = avg_reg, colour = region, shape = region)) + 
  theme(strip.text = element_text(
    color = "white"),
        strip.background = element_rect(
    fill = "grey20",
    color = "white"))

ggplotly(plotfacet)
cases_year1 <- cases_year %>% 
  group_by(year, region) %>% 
  mutate(avg_reg = mean(measles_incidence_rate_per_1000000_total_population)) 

hexgg <-  ggplot(cases_year1, aes(x =  year, y = avg_reg)) +
    stat_bin_hex(aes(fill = after_stat(density), colour = after_stat(density)), 
                 bins = 10,
                 linewidth = 1) +
    scale_fill_viridis_c(option = "B") +
    scale_color_viridis_c(option = "B", guide = "none") +
    labs(x = "Year", y = "Avg. incidence rate pr. region", fill = "",
         colour = "") +
    theme_minimal()
hexgg

plot_gg(hexgg, multicore = TRUE, windowsize = c(800, 800))

render_movie("silly.gif")
knitr::include_graphics('C:/Users/au483794/OneDrive - Aarhus universitet/Documents/Fridayproject/Peter/silly.gif')

cases_month %>%
  mutate(region = as.factor(region)) %>% 
  ggplot(aes(x = measles_suspect, y = measles_lab_confirmed, color = region)) +
  geom_point() + 
  facet_wrap(~region) + 
  theme_bw() + 
  labs(x = "Suspected cases of measles", y = "Laboratory confirmed cases of measles", color = "Region")

  cases_month %>%
  mutate(region = as.factor(region)) %>% 
  filter(region == "European Region") %>% 
  ggplot(aes(x = measles_suspect, y = measles_lab_confirmed, color = region)) +
  geom_point() + 
  facet_wrap(~region) + 
  theme_bw() + 
  labs(x = "Suspected cases of measles", y = "Laboratory confirmed cases of measles", color = "Region") + 
  coord_fixed()

library(dplyr)
library(ggplot2)
library(gganimate)
library(rnaturalearth)
library(sf)
library(magrittr)
# Load world map
world <- ne_countries(scale = "medium", returnclass = "sf") 
world <- world %>% 
  select(name, geometry)
#Calculating avg. incidente pr. 1,000,000 pr. country pr. year. 
cases_year_gif <- cases_year %>% 
  group_by(year, country) %>% 
  mutate(avg_country = mean(measles_incidence_rate_per_1000000_total_population))
# Joining data to map
map_data <- left_join(world,cases_year_gif, by = c("name" = "country"))
#Creating an interactive world-map showing incidence by country pr. year by color-coding. 
library(plotly)
fig <- plot_geo(cases_year_gif) %>%
  add_trace(
    z = ~avg_country,
    color = ~avg_country,
    colors = "Reds",
    text = ~country,
    locations = ~country,
    locationmode = "country names",
    frame = ~year
  ) %>%
  colorbar(title = "Incidence")
fig